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I ' Abstract 

^ I Using the parallel/orthogonal space method, we calculate the planar two-loop 

three-point diagram and two rotated reduced planar two-loop three-point diagrams. 
^ . Together with the crossed topology, these diagrams are the most complicated ones 

in the two-loop corrections necessary, for instance, for the decay of the boson. 
Instead of calculating particular decay processes, we present the new algorithm which 
allows one to calculate arbitrary NNLO calculations for massive planar two-loop 
vertex functions in the general mass case. All integration steps up to the last two 
ones are performed analytically and will be implemented under t^oo^s as part of the 
Mainz x^oop^-GiNaC project. The last two integrations are done numerically using 
methods like VEGAS and Divonne. Thresholds originating from Landau singularities 
are found and discussed in detail. In order to demonstrate the numeric stability of 
our methods we consider particular Feynman integrals which contribute to different 
physical processes. Our results can be generalized to the case of the crossed topology. 



1 Introduction 



Precision measurements at the LEP collider at CERN and other colliders like e.g. SLC at 
SLAC, TEVATRON at Fermilab and HERA at DESY have reached a precision which has 
exceeded all expectations. This is true especially for electron colliders [H |21 El El • At the 
moment the precision of measurements related to the parameters of the Standard Model 
of electroweak interactions reaches values up to 0(10"'^) 4J. At future colliders like the 
LHC or the ILC (including GigaZ), further improvements are expected jSl El El 13 IE] • 

Compared to this, theoretical predictions accomplish this precision only in very few 
cases. To check the validity of the Standard Model and to be able to draw conclusions about 
"new physics" , progress in theoretical methods and their application is necessary [3 El- The 
complexity of the calculations and the number of the graphs which have to be calculated 
within perturbation theory grows considerably order by order. At second order only a 
few observables are calculated Still there is no method available which allows for the 
(semi) automatic calculation of arbitrary processes at this order. Much work has been 
done on NNLO calculations |inilIIllI2llI3llI2liniliniinilIHl[ISI2nil2Ill23l2nil21 
ESIEHEZI- While for mixing QCD and electroweak 0{aas) NNLO vertex corrections a 
few calculations were done several years ago J71, for the evaluation of two- loop three- 
point diagrams in the general mass case only the methods proposed in Refs. (23 EOj 
have been used successfully in order to calculate electroweak NNLO corrections With 
increasing energy as it will be used at the ILC, radiative corrections will become increasingly 
important. Therefore, there is still need for independent methods to calculate general 
massive two-loop vertex diagrams. 

In this paper we present a new algorithm which allows one to calculate NNLO-correc- 
tions for general massive planar two-loop vertex functions as they arise in the Standard 
Model of electroweak interactions using Feynman gauge, e.g. for the effective weak mixing 
angle sin^ E21 occuring in Z° and for other processes like Z* tt. Even 

though the calculation of particular physical processes is not subject of this paper, the 
new algorithm enables one to perform such calculations. For our algorithm we use the 
parallel/orthogonal space method [33J which allows one to separate effects coming from 
inner momenta from those of momenta of the outer particles of the process. 

The benefit of the parallel/orthogonal space method is the fact that the calculation 
stays close to the physical process. The introduction of Feynman parameters and the ap- 
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plication of the Wick rotation which estrange the calculation from the physical process 
for other methods are not necessary for the parallel/orthogonal space method. Instead, 
a decay process can be calculated in the rest frame of the decaying particle. The intro- 
duction of Gram determinants which might cause artificial divergences is not necessary. 
The analytical integration leads to still rather simple basic functions like logarithms and 
dilogarithms. For the remaining numerical integration the integrand can be analyzed at 
physical thresholds. Landau singularities are mirrored directly onto the parameters of the 
integrands. Finally, the method is totally independent from other methods and therefore 
allows for an independent check. 

Up to now, the parallel/orthogonal space method was used successfully to evaluate the 
general massive scalar integrals in the case of the planar topology [34j and in the case of 
the crossed topology jSSj- However, it was not possible to evaluate two- loop three-point 
functions which contain loop momenta in the numerator. The method proposed in this 
paper allows one to calculate also tensor integrals for the planar vertex topologies. 

We present a tensor reduction which reduces the planar and rotated reduced planar 
two-loop three-point topologies containing an arbitrary tensor structure to a set of mas- 
ter integrals with strongly restricted numerator structure . As for most of the master 
integrals which remain after tensor reduction there are already established methods avail- 
able 1211 1221 12ni I2Z| , we do not give results in these cases. The same is valid for particular 
mass cases [TUl E7| EHl and for the UV-singular parts [Sn]- Efficient methods for several 
cases including few masses have been developed in Refs. [121 1201 121] • We also do not dwell 
on the calculation of massless diagrams for which other methods are applied successfully 
(see e.g. |^). As IR divergences occur in this case, there is no systematic method available 
within the parallel/orthogonal space method because for each special case the subtraction 
has to be found seperately "by hand" |H1 02] • Instead, we concentrate on the most com- 
plicated part, namely the calculation of the UV- finite parts of those master integrals which 
keep the planar or rotated reduced planar topology after tensor reduction. 

Two-loop integration in four-dimensional space-time needs eight integration steps. Us- 
ing the parallel/orthogonal space method, we can do six of them analytically using different 
techniques. The last two integration steps are done numerically. Because most of the inte- 
grations can be performed analytically, the Landau singularities of the Feynman diagrams 
can be related graphically to the integrands in the remaining two-dimensional integration 
region. Using different algorithms, we demonstrate the numerical stability. 
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The algorithm for the numerical calculation of the last two integrations is developed 
and its reliability is demonstrated. In addition, further exhaustive numerical tests are 
done. While for the numerical integration programs like the Monte Carlo integration 
routine VEGAS 03| ^^e program Divonne of the library CUBA |45^ are used, the 
implementation of the algorithms for the semi-automatical calculation in ;^[oops ITTl HH] 
is part of the work described in this paper. In the ;^[oops-GiNaC project developed by the 
members of the ThEP working group in Mainz |l^ l5Uj. work on the automatic generation of 
Feynman diagrams and their evaluation with analytical and numerical means is in progress. 
For this purpose the language for algebraic calculation GiNaC was developed While 
most of the master two-loop integrals not considered in this paper are already implemented 
in ^iCoops, the missing planar topologies will be added in the course of the work introduced 
here. The complete set of algorithms described in this paper is implemented under j([oops- 
GiNaC. However, the implementation of parts of the algorithms for known topologies as 
mentioned in Sec. 13.81 still has to be done. Because in this spirit the implementation is not 
yet finished, the calculation of physical processes is not subject of this paper but will be 
presented in a future publication. 

The paper is organized as follows. In Sec. 2 we introduce the main tools of the analytical 
calculation associated with the parallel/orthogonal space method. The tensor reduction 
procedure is explained in detail in Sec. 3 for the non-reduced planar two-loop three-point 
diagram, and modifications in case of the rotated reduced planar diagrams are mentioned 
and explained. The tensor reduction leads to master integrals which are integrated ana- 
lytically in Sec. 4 up to the last two integrations. In Sec. 5 we deal with the analysis of 
Landau singularities and thresholds. Sec. 6 is devoted to the numerical integration. In 
both sections we present examples to show the reliability of our procedure. In Sec. 7 we 
give our conclusions. In Appendix A we deal with the integral basis. 

2 Tools for the calculation 

In this section we provide the reader with the tools necessary for the calculation of mas- 
sive two-loop tensor vertex integrals. Most of the tools were already introduced in the 
literature [^21 OH 1121 1^21 1^21 IHl]- Therefore, we can be brief in presenting these tools. 
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2.1 The parallel/orthogonal space method 

The integrals we have to calculate are determined by the two momenta of the produced 
particles, pi and p2- Because the integrals are expressed in terms of covariant quantities, 
we are free to decompose any loop momentum k into two covariant vectors, k^ = k'^^ + kj_, 
where fcy has components in the parallel space which is the linear span of the external 
momenta pi, while k± is the orthogonal complement with k^pi^^ = with components 
in the orthogonal space |1U1 ITTj . For three-point functions we can consider the process 
in the rest frame of the decaying particle. In this frame the two emerging particles are 
produced back-to-back, defining the 2;-axis.^ In this frame the representation 

pi = (Ei; g„ 0, 0), p2 = {E2; -g., 0, 0), p = p, + p^ = (E; 0, 0, 0) (1) 

can be used where p is the momentum of the decaying particle with E = Ei + E2 and 
= E"^. Accordingly, the loop momenta k and / necessary for the description of the 
two-loop integral are parametrized by [SSI 

k = {ko;ki,k±), I = (WJiJi)- (2) 

The two-dimensional vectors k± and l± are represented in polar coordinates, using the 
squared absolute values s = fc^. = ^1 t = = and two angles, the angle a of k± 
with the X-axis and the relative angle 7. We can write 

k±-l± = k±l±z = y/siz, z = cos 7 = cos(fc_L, (3) 

The choice of parallel and orthogonal subspaces itself is Lorentz invariant, so that the 
calculation can still be done in any Lorentz frame. The benefit of the parallel/orthogonal 
space method (P/O-space method) is that the contributions belonging to the orthogonal 
space can be integrated out and we are left with the contributions in parallel space only. 
The integration measure is accordingly written as 

d k d I = TT dko dlo dki dli / dsdt l —j^=^= (4) 
J J -00 Jo J-i \/\ — 

where the trivial integration over a has already been performed. 

^The z-component will be written as the second component of the four- vector in the following. 
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2.2 The linearization 

Using Feynman gauge, the denominators of integrals occurring in two-loop vertex calcula- 
tions contain up to six propagator factors. A typical factor is given by 

Pi = {k + pif -ml + ir] (5) 

where rj > In the P/O-space representation this propagator factor reads 

Pi = (A;o + E^f - (fci + q,f -kl-ml + tr]. (6) 

If we replace /cq = ± ^i, we obtain 

Pi = (k', + ± 2k,{k', + E,T Qz) - ql - kl - mj + it] (7) 

where fci no longer appears quadratically. This replacement, also known as linearization, 
is allowed because the integral for ko ranges from — oo to +oo. However, because of the 
occurrence of the mixing propagator factor 

Ps = (A; + ly -ml + IT], (8) 

the signs for the linearizations in k and / are coupled. In applications of this linearization 
we use the sign which is the most appropriate for our aims. The benefit of the linearization 
is that the integrals over ki and h can then be calculated by using the residue theorem. 

2.3 The integration over z 

The quantity z occurs only in the just mentioned propagator factor P3. After linearization 
this propagator factor can be written as P3 = A + Bz + irj where 

A = A{ki,li) = {k'^^ + Q''±2{ki + li)% + Q-{s + t)-ml, B = -2^t. (9) 

P3 occurs only once in the denominator (if at all). The resulting integral 

dz \ 71 TT 



-1 VI - A + Bz + iri ^{A + ir]Y-B^ d{A + i7]Y - Ast 



(10) 

can be written in this closed form only if we make a convention different from the usual 
one |TUj. We stipulate that the cut of the square root is located on the positive real axis 



instead of the negative one. For this reason we have added the star in d{A + iiqY — P^ 



which reminds one of the different cut. The occurrence of the modified square root has 
consequences for the subsequent integrations because in closing paths for the residue the- 
orem we have to avoid crossing this cut. In order to find conditions for this we spht the 
square root up into a product of two square roots \ri¥±^~+lr]. If one of the radicands is 
equal to a positive real number x, we are definitely crossing the cut. Written in terms of 
ki and li, we obtain 

^'^'^ ~ ■ ^^^^ 

Note that the two undetermined signs are not correlated. While the sign between y/s and 
Vi depends on which of the two square roots is taken, the global sign of the right hand side 
is due to the sign of the linearization. Considering only the linearizations ko — > k'^ + ki, 
lo ^ I'o + h, the signs of the real and imaginary part of the cut still depend on the sign of 
(/cq + I'o). If (/cq + I'o) > 0, the cut is located in the lower complex half plane. The real part 
starts (for x — 0) at some finite value and runs to +00, independent of which square root 
we selected. In order to avoid the cut, we therefore have to close the integration path for 
ki resp. li in the upper complex half plane while for (/cq + /q) < we have to close it in the 
lower half plane. The situation is opposite for the linearization with a minus sign. 

3 Tensor reduction 

After having introduced the main tools for the calculation, we can start with the calculation 
itself. The starting point is the tensor integral 

H^i)-/^°/?H/i)''^(A;xM^4, ,4, 

aoaia2bobib2C — J P1P2P3P4P5PQ anat {IZ,) 

where 

Pi = {k+ PiY - 1711 + if] Pi = {I - Pif - ml + ir] 
P2 = {k-p2f -ml + iri P5 = (l+P2f -ml + ir] 
P3^{k + If -ml + ir] PQ^f-ml + ir]. (13) 

If one is doing Standard Model calculations using Feynman gauge, the powers are restricted 
by < 0,0, a-i, 2a2 < 3 and < h^, bi, 262, c < 4, as well as by < ao + Oi + 2a2 + c < 3 and 
i < bo + bi + 262 + c < 4. Nevertheless, the algorithm introduced in this paper works for 



7 




Pi 



p 



■> 



P2 



Figure 1: Momentum flow convention for the planar two- loop three-point function 

general powers. In Fig. ^ we show the momentum flow for the diagram in order to define 
the momenta. For the indices a^, 6j, c we allow non- negative integer values. In performing 
the tensor reduction the integrals are simplified to integrals with simpler numerator and/or 
denominator structure. 

For one-loop integrals the numerator can always be removed by reduction procedures 
[S3]. For two-loop integrals this need not be the case. In general there are not enough 
propagator factors to cancel all components of the numerator that occur. For the genuine 
planar two-loop three-point function all numerator factors related to the second loop mo- 
mentum / can be cancelled. This is not the case for the reduced topologies, and it is not 
the case for the first loop momentum k. The aim of the tensor reduction in general is to 
reduce the numerator as far as possible so that the integrations can be performed similar 
to what is done for a trivial numerator The reduction procedure will be explained in 
the following. The representation of the procedure by diagrams will show the topologies 
only. In this spirit plus signs like the one occuring in Eq. (fTU|) have to be understood as 
sums over diagrams with the same topology but different factors and signs. 

3.1 Cancellation of the mixed contribution 

The first step in the cancellation procedure consist in cancelling powers of the mixed factor 
{k±l±z). This factor occurs also in the propagator factor P3. We can write 



(14) 
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where 

N^^kl-kl-k]_ + ll-ll-fj_ + 2kolo - 2^iZi -ml + irj. 
Used iteratively, we obtain 



{kj±zy 



^ ^ i=0 ^ 



(15) 



(16) 



This iterative formula is necessary because the propagator factor P3 occurs only once. As 
applied to the integral, we obtain 



{kJ^zfd^kdH _ 1 r N^d^kdH 
P1P2P3P4P5P6 ~ P1P2P3P4P5P6 



c~i—l 



y^J_ f Nljkj^zj 
^ 2^+1 J P,P2P^P5Pe 



-d^kdn 



(17) 



where we have skipped all the other numerator factors for convenience. The first part is 
the same planar integral with a different numerator structure. However, the second part 
no longer contains the mixing propagator factor P3. Instead, the diagram factorizes. In 
addition, we can conclude that this vanishes if c — z — 1 is odd. The reason is that the 
integration over z is given by 



+1 z^-''-^dz 



-1 y/T^ 

which vanishes if the integrand is odd. Diagramatically we can write 



(18) 




X {kj,±zy 




+ 




(19) 



After having cancelled the mixing part it is irrelevant with which loop momentum we 
continue. However, we always start with the orthogonal space components. 



3.2 Cancellation of 

The factor k\ — s occurs in the denominator factor Pi. We can write 

k'i = Ni- Pi, Ni^kl-kl + 2koEi - 2kiq, + -ql-ml + irj. 



The iterative formula 



a2-l 

(kl)"' = iVf - J2 N{{klY^-'-^Pi 

i=0 



(20) 



(21) 
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can be used for the integral to obtain 



-P1-P2-P3-P4-P5-P6 J P1P2P3P4P5P6 
Schematically this reduction reads 

4 



02 — 1 

E 

i=0 



^-^^ d'kdH. 

-P2-P3-P4-P5-P6 




+ 

5 5 

In the following only the diagrammatic reduction will be shown. 





(22) 



(23) 



3.3 Replacement of ki 

In anticipation of the linearization we replace the numerator factors ki by ^{{ko^ ki) — ko) , 



usmg 



i=0 



(24) 



While the additional factors /cq^ * will be cancelled together with kQ° in the next step, the 
subsequent linearization will replace {ko T ki) by k'^. 



3.4 Cancellation of ko 

Pi and P2 are the only propagator factors that contain only k and not We can combine 
these both to obtain 



2Eko ^ N2 + Pi - P2, N2^pl-pl-ml + ml 
An iterative formula can be constructed as before, the reduction reads 





+ 




+ 



(25) 
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3.5 Cancellation of / 



_L 



In using 



li = Ne- Pe, Ne = l'o-li-mi + ir] 
we can cancel the factors ending up with the reduction scheme 



+ 



(27) 






(28) 



3.6 Cancellation of li 

Different from the situation for the replacement of /ci, we have three propagator factors 
which do not contain k. In this case the factors l\ can be cancelled completely. We use 



to reduce powers of h according to the scheme 




6 Xll' 




+ 




+ 



(29) 




3.7 Cancellation of Iq 

Finally, we can cancel the factors Iq as well, using 

2Elo = iVs + P5 - P4, N5=pl-pl-ml + ml 

The reduction scheme reads 

4 




6 X 



^0° 




+ 




+ 



(31) 
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3.8 Summary of the reduction procedure 

After performing the reduction procedure, we can collect all the steps into a single one. All 
diagrams which are still of the (original) planar topology after reduction procedure have a 
numerator factor {k^ =F^i) in different powers. The maximal power, i.e. the maximal value 
for a, is given by ai+2a2+2c. The numerator factors of the diagrams with reduced topology 
are changed but not easier. Nevertheless, because of the different topology, reduction 
procedures which are known in the literature or which will be explained in the following 
can be applied. Schematically we write 




Denoting the topologies by T^,T^, . . . , according to the cancelled propagator factor, 
we can cite the following references: 



• Two-loop three-point functions with two-point subloop, as they are given in our case by 
the topologies and T^, can be calculated by using dispersion relations PH 123 12^] ■ 

• The factorizing topology is calculated as a product of one-loop topologies. This 
contributions may be calculated using ?([oops |50j or similar packages like SANC [321 • 

• The effective two-loop two-point topology can be calculated using ?([oops |77j . 

• What is left are the rotated reduced planar topologies and T^. They can be reduced 
in a similar manner as the original planar topologies.^ In the following section we will 
deal with these topologies and the necessary modifications of the reduction procedure. 

^Because of historical reasons, these topologies are called rotated reduced planar topologies. In the 
original reduced planar topology the line crossing the triangle is going horizontally from left to right, 
connecting the left vertex with the middle of the vertical line on the right. 
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3.9 Modifications for the rotated reduced planar topologies 

For the rotated reduced planar topologies we assume the same numerator factor as for the 
non-reduced one, ko^'k^' {klYH^^nl' {llf^{kJ±zY. The procedure is quite the same as m 
the case discussed before. But because one of the denominator factors P4 and P5 is absent, 
it is not possible to cancel the factors of li as in the previous case. Instead, we anticipate 
the linearization as in the case of ki (and with the same sign) in writing 

h^T ((/o T h) - /o) ^ ll' = iTlf' j:i-lo)''-' (^'] ik T hy. (34) 

Again, while the additional factors Iq^~^ will be cancelled together with 1^° in the last step, 
the subsequent hnearization will replace {Iq ^ h) by I'q. However, this last step makes the 
difference between the two topologies. 

For the topology the denominator factor P4 is absent. However, we can obtain 

2{E2 ± q,)k ^K-Pe + K = ±2(/o T h)qz - P2 - ^6 + ml (35) 

At this point the different signs for the linearization enter the game. While we have 
designed the term {Iq ^ li) in the way that it can be combined with the factors from the 
cancellation of li, in order to obtain Iq we have to divide by (£'2 ± Qz)- However, if the 
second produced particle is massless, p| = 0, this factor will become zero. In order to 
include also the massless case, we will prefer the linearization li — l[ + Iq for this topology, 
i.e. the upper sign. The whole reduction procedure for the case of the topology reads 
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All the diagrams on the right hand side except for the first one with maximal powers 
Oi + 2a2 + 2c for a and 60 + + 262 + 2c for are again members of the different simpler 
topology classes which were mentioned in the previous paragraph. 

Looking at the topology T^, the propagator factor P5 is absent. In this case we can 
obtain a reduction formula from combining P4 and Pq to obtain 

2(^1 T q.)lo = K-P^ + Pe, K = t2(/o T h)qz +pl-ml + ml (37) 

In this case the linearization li = l[ — Iq (lower sign) is more appropriate. The result of 
the reduction reads 




4 Integration of the master integrals 

After having reduced the numerator of the integrand to factors {ko =1= ki)°' in case of the 
planar two-loop topology and (/^o =F ^i)"(^o T^i)'^ in case of the rotated reduced planar two- 
loop topologies, we can start to integrate this set of master integrals T^, 7^^^. Before we do 
so, we have to consider the occurrence of ultraviolet (UV) divergences and their treatment 
in terms of appropriate subtractions. As for the treatment of infrared (IR) divergences we 
refer to Ref. PT] for the scalar case. In the more general case, methods taken from Ref. jSZI 
can be applied. However, in this paper we consider only IR-finite examples. 
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4.1 UV-divergences and subtraction procedure 

The UV-divergences can be subdivided into three classes: divergences with respect to the 
two loop momenta and global divergences. It is quite obvious that if we have dk propagator 
factors depending only on k, di propagator factors depending only /, d^i propagator factors 
depending both on k and /, and a power A;"*^/"' in the numerator, the corresponding degrees 
of divergence in D space-time dimensions are given by 

iVk = 2(4 + dki) -Uk- D, 
uji = 2{di + dki) -ni - D, 

uj = 2{dk + di + dki) - (nk + rii) -2D. (39) 

The two-loop integral is divergent in D dimensions, if at least one of the degrees u^, uji, or u 
is zero or negative. If the degree of divergence is zero, this divergence is called logarithmic 
divergence. 

In order to calculate the integral, we first have to regularize it. We use dimensional 
regularization and write D = A — 2e. After that, integrals can be split off into a divergent 
and a convergent part. While the convergent part can be calculated for e = 0, i.e. for D = 4 
space-time dimensions, the integrand of the divergent part is simpler but has the same UV- 
behaviour. In an appropriate subtraction procedure, therefore, we subtract and add an 
integrand which is simple enough to be integrated analytically (at least the singular part) 
but has the same singularities as the integrand itself in order to cancel the singularities. 
This subtraction procedure is found and will be formulated multiplicatively |32] ■ The 
subtracted integrand contains the subtraction factors 



with 



d 



1=1 i=l 
di di 

A = n ((^ + Pl.^y - + , A,r- = n ((^ + ^rP^f - ml,r + ^V) (41) 

1=1 i=l 

where Pk,i,Pi,i e {±pi,±p2}, rrik^i G {mi,m2}, mi^i G {m^^m^}, and mk,i,r,mi,i,r are sub- 
traction masses which have to be introduced artificially.^ The degrees of divergence for this 



^In the following wc will simplify the indices according to the special cases. 
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subtracted integrand increase by jk and ji, resp. and thus may lead to integrals without sin- 
gularities for arbitrary values for Kr, Xr, 'rrtk,i,r-, and mi^i^^. Finally, the limit k^, — ^ can 
be performed, leaving the subtracted integral convergent and enabling the semi-analytical 
calculation. Since the subtraction terms always contain a two-point subloop, the singular 
part can be calculated analytically. For the non-singular part one can use e = 
In the end the sum of subtracted integral and subtraction terms have to be independent 
of the subtraction masses. In the following we will deal with the subtracted integrals only. 



4.2 Subtraction for the planar topology 

The reduction procedure explained in the previous section leads to the integrals 

The degrees of divergence can be calculated in four space-time dimensions, 

tUfc = 2 — a, LUi = 4, uj = A — a. (43) 

There are no sub divergences in /, i.e. divergences caused by the integration over the loop 
momentum /. However, for a > 2 we have to subtract subdivergences in k. We choose 
jk = a — 1 and multiply the integrand in Eq. ()42|) by 



where Pi = {k + — ml + irj and P2 = {k — ~ ^2 + ^^7 given in Eq. (fT!?|l . The 
modified propagator factors read 

Pi,i = ik + KiPif -mli + ir], P2,i = {k - KiP2y - ml , + ir]. (45) 



Note that each of the subtraction factors in Eq. (j44j) improves the degree of divergence at 
least by 1 because by power counting the denominator is given by fc^ and terms of lower 
order in k while the numerator starts with k^. The masses mi^i,m2,i of the subtraction 
can be chosen arbitrarily. However, in order to avoid the introduction of spurious IR- 
singularities and squared propagator factors in the denominator, the subtraction masses 
are usually chosen to be non-zero and different from the physical ones and from each 
other 02] ■ Furthermore, in order not to introduce new thresholds, they should be larger 
than the masses of the decaying particles. On the other hand, they should not be too large 
because of the possible occurence of numerical instabilities. 
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4.3 Integration for the planar topology 

After the linearization fco = k'^ + ki, in the P/O-space representation the propagator factors 
are given by 

2{k'^ + Ei- q,)ki + {k'^ + Elf -ql-s-ml + iiq 
2{k'o -E2- q,)ki + {k'o - E2f -ql- s-ml + iiq 
2 (fco - K.i{Ei - q^)) ki + (/cg + KiEiY - K^g^ - s- + irj 
2 (fc; - /€i(E2 + q,)) ki + (A;; - i^iE^f - K^ql - s - + i7]. (46) 

In case of the (non-reduced) planar diagram where both linearizations are possible, we 
select the linearizations k^ = k'^ + ki and = I'q + h for convenience. The propagator 
factors which are linear expressions in ki can be written in closed form expressions, 

Pi = 6iki + - s + ir], Pi,j = 6i^iki + C.i,i - s + ir] 

P2 = 02ki +^2-s + ir], P2,* = 02,iki + ^2,^ -s + ir]. (47) 

Because the parameters can take arbitrary values, the denominator of IC^^''^ in Eq. (j4(J|) 
together with the factors Pi and P2 from the integrand consists of different linear factors. 

Integrating over ki, the integration path [—00, +00] can be closed in the upper or lower 
complex half plane, and Cauchy's theorem can be applied. The decision for one of these 
paths depends on the sign of the linearization and on the sign of (/cq + I'q) (cf. Eq. 
and its discussion). The poles are given by the (different) zeros of the propagator factors 
present in the calculation, the residues are given by 
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(48) 



The condition Pi 
1 

^2 Pi=o ( 



results in ki 
Oi 



[s — ^1 — irj)/ 61 =: A;i(i). One obtains 
A/? „ _ 02i 



dl)is - S2,l - ii]) 



S2,l - IV 



S2,l 



O1C2 



02 — dl 



(49) 
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In the same way P2 — results in ki — {s — ^2 — iv)/^2 —'■ ^1(2)- One obtains 



1 



Ol^2 — 02^1 



(50) 



P2=0 (^1 - ^2)(S - Si,2 - 8-31^2-111 ' -6*2 

Note that si,2 = S2,i- For the other two types of residues, however, a simphfication appears 
after inserting the poles. Because of the fact that for the semi-analytical calculation of the 
convergent part of the integral we use the limit Ki ^ for all i, we obtain 



2A 



2ko-.9o for i = l,2,...,jfe. 



(51) 



In this limit, and on the pole Pij = given by ki = {s — ^ij — if]) / Oq =: k^ij), the different 
propagator factors read 



P2,i 
Pi 
P2 



P,,.=o = ^((^i-^o)(«-ir?) + (ei^o-^i6j)), 

p^^^^O = ^((^2-^0)(«-i^) + (6^0-^26,,))- 



Obviously, the dependence on s appears only in the numerator. We can write 

jk / 



-ft 

°2j i=l V Pl,iP2,i, 



and similarly 
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.Jk 



n 1 



The different poles are given by 



P2,,=0 



^i(i) 
ki{i,i) = 



P1P2 

Pl,iP2,i 

s-i2-iv 



2jk-2 

0=0 



2jk-2 
a=0 



(52) 



(53) 



(54) 



01 

s - ^i,i - iv 



ki(2) = 



a ■~2{2,i) a 

For < (i = 0, 1, 2) the corresponding pole can be found in the upper complex half 
plane. For the hnearization ko — k'Q-\- ki and k'^ + l'^ > 0, the path has to be closed in the 
upper half plane. In this case we obtain a non- vanishing residue. For /cq + ^0 < 0, however, 
the residue occurs only in case of > and has the opposite sign. Using 



O2 



IT] 



for i = 1,2, . . .,jk- 



(55) 



B'^ = 2m [e{k', + Qe{-On) - 0{-{k'o + QWiOn)] , 0,1,2 



(56) 
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where 6{x) is the step function, we can integrate over ki to obtain 

:^(n(-a;))^'-) 

2 / Kk\rk 3k 2ife-2 

= E ^^^fih^n). ^l) + E E KK.^aS'fihin,), h) \ (57) 



where f{ki, h) is given e.g. by the inverse of ^ {A{ki, li) + irjy — 4st (cf. Eq. flTUj)). 

For the integration over li we have to consider the propagator factors P4, P5, and Pq. 
Because the degree of divergence ui is positive, no subtractions are needed. After the 
hnearization /o = ^0 + ^1 '^^^ write 

P4 = 2{l'o-Ei + q,)k + {l'Q-Ei)^ -ql~t-ml + tr], 

P5 = 2(/[, + P2 + g.)/i + (/; + ^2)'-g'-t-m2 + 2r/, 

Pg = 2/[,/i + - t - + (58) 

which again can be cast into the closed form expression 

p. = 0./^ j^^.-t + irj, i = 4, 5, 6. (59) 

The procedure which was explained for the integration over ki works accordingly. The 
poles of 

^ fn(M+^i-^+^^)j (60) 



P4P5P6 

are found at 



t — 1pm — if] 

h{m) = 7 form = 4, 5,6. (61) 

From Eq. (jllll we read off that the cut is located in the same complex half plane as in 
case of the integration over ki. Therefore, the occurrence and the sign of the residue is 
determined by 

= 2m [e{k', + l'o)0{-(pn.) - e{-{k'^ + Q)e{<Pm)] , m = 4, 5, 6. (62) 
Using Cauchy's theorem, we obtain 

/ + 00 rl]^ 
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E 

m=4 



6 

E 

m=4 



n WUm) + -t + 17]) 



j=4 



(63) 



k2 Ul 



Of course, the same is valid if ki is replaced by kun) or ki(n,i)- Usually, (pj and 0m are 
different for j ^ m. In this case we can write 



^(0i - 0m) - (0i^m - 0mV'j) " ^^7(0^ " 0m) 



- 0m) - ^j,m - «'?) 



with 



0jV'm - 0mV'i 



We perform a partial fraction decomposition to obtain 



m=4 j=m+l ''mj '"'/ 



i3' A/"' ■ 

f{ki,l 



l(m), 



(64) 



(65) 



(66) 



-oo P4P5P6 

where A/'^.-,- are the corresponding coefficients including 0^. A special situation occurs 
if = 0. Because of Qz = Ei, in this case we have 04 = 06 and the corresponding 
denominator factor contributes in the form ("04 — 0'6)06 only. However, if P2 = 0, no special 
case has to be taken.'* 

The integrations over ki and li and the integration over z can be combined in the closed 
form expression for the convergent part V° of the master integral in case of Hi 0, 



jko - k,Y 



-JC^'^^d^kdH 



2 5 
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-^EE E 

n=l m=4 j=m+l 



00 A/"' ./7f 



t — t 



IT] 



+ E E 



4st 



(67) 



where 



A:/ 



-47r2 [^(fc; + /;)0(-e„)^^(-0„J + e{-{k'^ + /^,))^(^J^(0m)] (6J 



and An,m,An^i^m are the coefficients A in Eq. Q with ki replaced by fci(n), fci(n,j) and li 
replaced by i.e. An,m = A{ki(^n),h{m)) and A„,i,m = ^(fci(n,i), ^i{m))- 
^For the other hnearization lo = I'q ^ h the situation is the opposite. 
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4.4 Subtraction for the rotated reduced planar topologies 

For the rotated reduced planar topologies we are left with the integrals 

where in case of the topology (x = 5) the upper sign, in case of the topology (x = 4) 
the lower sign is valid. The degrees of divergence can be calculated to be 

ujk = 2-a, uJi = 2-f3, u = 2-a-p. (70) 

For (3 > 2 the integration over I becomes divergent and a subtraction has to be performed 
in addition to the subtraction in k. The subtraction factor is given by 

C^^'^ = fl (l - p-] , x = A,5 (71) 

j=l \ / 

where in addition to P4 = (/ — pi)^ — m\ + irj and P5 = (Z + — + ir] we have 

^4,i = {l- \jPi? - rri\^ + ^r?, = (^ + hV^f - m\. + ir]. (72) 

The same limitations as they were mentioned at the end of Sec. 14.21 for the subtraction 
masses mi j and m2,i apply for the masses m^^j and ms j as well. 

4.5 Integration for the rotated reduced planar topologies 

The integration over k\ can be done in the same way as for the original (non-reduced) planar 
topology. However, special care has to be taken concerning the different linearization in 
case of the topology T^. For the topology we have = — 2/cq and the step function 
term has to be replaced by 

= 2m [e{-ik', + Q)e{-e^) - e{k', + Qe^er^)] . (73) 

For the integration over li we have to perform a subtraction. After the linearization 
/o = ± ki the propagator factors in the P / 0-space representation are given by 

P5 = 2{l'^ + {E2 + q,))h + {l'Q + E2f -ql-t-ml + ir] 
P5j = 2{l'^ + \j{E2 + q,))h + {l'^ + \jE2f -\]ql-t-ml + ir] 
Pe = 2/oZi + Iq ~t-ml + ir] for the topology T4, (74) 
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P4 = 



-2{1'q - {El + q,))h + {I'o - Elf -ql-t-ml + ir) 
-2(/[, - X,{E2 + q,))h + {I'o - A,^2)' - Xhl - t - 
-2/0/1 + I0 - t-ml + iri 



for the topology T^. (75) 



These propagator factors can be written as 



Px = 0x^1 + iIjx -t + ir], Px,j = (f)x,jh + i^x,3 -t + ir] x = 4, 5 

Pe = 06^1 + ^l^6-t + ir]. (76) 



In the hmit Aj — > 0, which will be considered in the following, we see that 



(77) 



In case of the linearization /q = I'^ + h the integration path [—00, +00] over li can be closed 
in the upper complex half plane for (/cq + /q) > and in the lower complex half plane for 
(/cq + /q) < 0. For the other linearization = k'^ — ki the situation is just the opposite. 
Because the poles arc given by 

t-ip^-ir] , t- ip^j - IT] , t-ipQ-ir] 



l{x) 



the integration ranges in t are constrained by 

^^ ^{2m[e{-{k', + l'^))e{-(t>^)-e{k'^ + l'^)e{ct>^)] for x = 4 (topology r^) 
1 2m [e{k'^ + Qe{-^^) - e{-{k', + /[,))^(0™)] for x = 5 (topology T^) 

In order to use Cauchy's theorem, we calculate the residues 
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For the first residue we obtain 
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(80) 
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Px=o (06 - (t>x) {t - te^x -iv) t- h,x - irj ' 
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''6,2; 



06'0a: " 0a:V'6 
06 - 01 



(81) 



For the second one we first calculate 



^6 - ^x,j, 



1 



Px,j=0 



(82) 



Therefore, the denominator of this residuum is independent of t while the numerator is a 
power series up to the power P'-~^. We define the coefficients A/^j^a by 



6 i=l 



(83) 



The last residue appears to be a combination of a pole in t and a power series. However, 
we can separate these two parts by adding and subtracting an appropriate term. This term 
is given by 



1 



tx.a — 



(l)xA - Mx 



Ps=Q ~ ^e) {t - t^^Q -ir)) ' t- t^fi -IT]' ^ 

If we subtract this term from the result for the last residue, we obtain 



(84) 



1 



p. 



^ = M^-t^ + 0{Pl)\-^ = -t^ + 0{P.). (85) 



a; i=l 



^ p 



P 



The difference no longer has a factor P^ in the denominator. Instead, we obtain again a 
power series up to the power We define the coefficients J^ij^a by 



1 A / Px 



X i=l 



P6=0 



(86) 



Having calculated the residues, we can perform the integration over li to obtain 
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The integrations over ki and /i and the integration over z can be combined in the 
closed form expression for the convergent part V^~^ of the master integral T^~^ in case of 
Ki, Xi 0. One has 



X 



2 r Qkl j^kj^l 



E E 

n=l m=x,6 I (S - Si,2 - iVj) {t - t^fi - 17]) ^ {An,m + irjf - 4:St 
jk 'ijk-'i nkl Kfk a"- Kfl 



i=l 



h jl-^ Ukl KfkKfl fb 



j=i 6=0 (s - si,2 - iri) ^ iAn,m,j + iriY - Ast 

jk jl 2jk-2ji-l nkl \fk \fl fb ~\ 

+EE E ^'j<>> [ (^gg) 

i=lj=l a=0 6=0 y (^n,i,mj + ^^7)^ " 4st J 

For the topology (a; = 5), B^''^ is given by Eq. (jUH)) . i.e. it is the same as for the 
non- reduced planar topology. For the topology (x = 4), however, we obtain 

Bt[^ = B'A = -4vr^ [^(-(^0 + lom-On)e{-<P,n) + ^(A:(, + Qe{e^)e{<Pm)] . (89) 
The coefficients in the square roots read A^^m = A{ki(n),l^rn)), An,i,m = A{k^n,i),h(m)), 

An,m,j = A(ki(^n^,li(^m,j))y SHcl An^i^m,j = ^(^l(n,i) ; ^(m,j) ) • 

4.6 Integral basis for the orthogonal space quadrature 

The integrals over s and t resulting for both the non-reduced and the rotated reduced 
planar topology are called basic integrals. They are of the types 

/""o dsdt 

T{rs,rt,ro,So,to) 







[s - So - if]) {t~to- ir]) ^ {tq - r^s - nty - Ast 
s^ds dt 











(t - to - iv) V (^0 - rsS - Ttty - Ast 
t^ds dt 

{s - So- ir]) ^ (ro - r^s - nt)^ - Ast 



K,irs,ruro) = / , sH^dsdt ^^^^ 

^ [tq — r sS — Ttt)'^ — Ast 
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All these integrals are divergent except for the first one. However, after the subtraction 
they occur in sums in which the divergences vanish. While the integral jF(rs, r^, tq, Sq, to) 
is calculated in Refs. IHEl IMl EH] , the other integrals are new. The results can be found 
in Appendix 1X1 in terms of the parameters r^, r^, tq, sq, and to- In this paragraph we are 
dealing only with the dependence on the parameters r^, r^, and tq and on the poles found 
in the residue integration. In case of the integral J-'^r^, rt, tq, Sq, to) we obtain 

To = [ko + /o) T 2 — + —7 (fco + ^o) -m^-iV =■ rn,m, 

Ts = 1T2^, n = 1^2^ (91) 

where the signs correspond to the two linearizations Iq = l'^ ± li. Note that r^, > and 
VsTt > 1 in regions where the numerical integration has to be done. For JF* we have to 
replace 6'„ by 9q and ^„ by ^n,i- For JF^ we have to replace 0^ by (j)Q and ipm by ipmj- For 
J-'^/s, finally, both replacements have to be performed. Note that up to the choice for the 
parameters so, to, and rt, J-"^ and are the same integrals. 

The imaginary parts in the rational factors in the numerator are used to separate real 
and imaginary parts according to the Sokhotsky-Plemely relations 

f(x)dx ^ f(x)dx . ^, . r, x , ,r.r.^ 
hm / ■'^ ' = P/ — \ S(x - xo)f(x)dx (92) 

V^+O J a X — Xq ± ITj Ja X — Xq J a 

where "P" indicates the principal value integral. 

4.7 Analytic behaviour of the modified square root 

We have to consider different cases in order to analyze the analytic behaviour of the square 
root occurring in the integrand [SH ESI El] • The equation 

R{s, t) = (ro - r,s - ntf - 4st = (93) 

parametrizes an ellipse which separates the positive and negative values of the radicand. 
Outside of the ellipse the radicand is positive, while for points inside the ellipse it takes 
negative values. The ellipse touches the axes at s = ro/r^ and t = ro/r^. For ro < 0, 
therefore, the ellipse is located in the third quadrant. But because the integration is 
performed in the first quadrant (s G [0,oo], t G [0,oo]), the square root will give real 
values and no imaginary part is produced. For ro > 0, however, the ellipse moves into the 
integration region and imaginary parts enter the calculations. 
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Figure 2: Regions in the (s, t)-plane with different real and imaginary parts for the radicand 

R{s, t) = (ro — TsS — rttY — Ast 

In order to evaluate the square root uniquely, we also consider the imaginary part of 
the radicand, caused by the term irj. It turns out that the sign of the imaginary part is 
determined by the sign of tq — r^s — r^t, the zeros of which are given by the straight line 
connecting the points s = tq/ts and t = ro/rt on the axes where the ellipse touches the 
axes. Above this line (as seen from the origin for ro > 0), the imaginary part is negative, 
below this line it is positive. The situation is shown in Fig. |21 

We distinguish the following cases: 

• outside the ellipse and below the line (Rei?(s,t) > 0, Im_R(s,t) > 0) we have 

^ (ro - TsS - nty - Ast = (ro - r^s - - Ast 

• outside the ellipse and above the line (Rei?(s,t) > 0, Im_R(s,t) < 0) we have 

\/ (ro - - nty - Ast = (ro - r^s - nty - Ast 

• inside the ellipse (Re R{s, t) < 0) we have 

\/ (ro - TsS - Ttty - Ast = i^jAst - (ro - r^s - Vtty 

5 Landau singularities and thresholds 

Multiloop integrals with many propagator factors are influenced by their singularity struc- 
ture. The pole structure can be analyzed by using Feynman parametrization [SHI 1^ . 
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Starting from 



F{pj, mt) = / fn ^^^') /fe' "^k) n ^n = ql-ml + iri (94) 



^i=l / n=l 



with £■ outer momenta pj, L loop momenta ki, N inner masses m„, and g„ a linear com- 
bination of Pj and A;;, we can apply the Feynman parametrization 



1 



AiA2---A 



N 



^ ^'io (Ai/li + A2A2 + . . . + A^A^)^ ■ ^ ^ 



By shifting the inner momenta ki to k[ by finite amounts (which are determined by the 
outer momenta), the denominator can be rewritten as 

N 

V=J2 >^nAn = + K{k[, k'^,..., k'^) (96) 
n=l 

where is independent of inner momenta and if is a quadratic form in the shifted momenta 
k[. It can then be shown 150] that the denominator vanishes, i.e. the integral has so- 
called Landau singularities, if 

A„(g,^ — m^) = for all n = 1, ... (first Landau equation) (97) 

N 

and ^ A„g„ = for each loop (second Landau equation). (98) 

n=l 





Figure 3: Sunset diagram with momenta (left hand side) and masses (right hand side) 



5.1 Two- and three-particle thresholds 

The Landau equations can be used to analyze a given Feynman diagram. For two-point 
functions, singularities occur if the squared outer momentum crosses thresholds or pseu- 
dothresholds. For a genuine sunset diagram as shown in Fig. 13 a two-loop two-point 



27 



function with three hnes of mass values mi, 1712 and connecting the incoming with the 
outgoing leg, these thresholds and pseudothresholds are given by 

= (mi ± m2 ± (99) 

where the three-particle threshold is located at = (mi + 1712 + ma)^. In general, the 
genuine threshold is the one with the highest value. If the squared momentum crosses this 
threshold, an imaginary part appears. According to the Cutkosky rules [HI], this imaginary 
part corresponds to the situation where the particles move onto their mass shell. In terms 
of Feynman diagrams, inner lines are cut into outer legs for physical particles. This can 
be visualized in Feynman diagrams by drawing a cutting line. 

In the diagrams we are working with we expect two- and three-particle thresholds. 
While the two-particle thresholds are related to the values of Sq = Sn,m and to = tn,m, all 
the three-particle thresholds are related to the values of tq = rn^m, 
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■Si,2 > ^ > [nil + m2Y , ''^1,5 > <^ > {mi + + m^Y , 

^4,5 > <^ > {m^ + m^Y , ^2,4 > <^ p^ > {m2 + + m^Y , 

^4,6 > <^ p^ > {m^ + mof , ri,6 > <^ p\ > {mi + m-i + mof , 

hfi > <^ p^ > {m^ + m^Y, r2,6 > ^ pj > {m2 + m^ + m^^f. (100) 

Anticipating examples to follow, for the subtraction terms some of the two-particle thresh- 
olds vanish while new three-particle thresholds appear. For the diagram shown on the right 
hand side of Fig. additional three-particle thresholds appear for p\ > {mnj + ms -|- m4)^ 
('^nj,4 > 0) and P2 > {rrinj + m^ + rn^Y i^n,j,5 > 0) while two-particle thresholds including 
the masses mi and m2 no longer occur. 

5.2 Landau singularities for the scalar topology 



As an example we consider the original planar topology in Eq. (P2j) with a = 0. Because 
we are only interested in a qualitative analysis, we use the following fictitious values for 
masses and outgoing momenta, 

mi = 420 GeV ms = 100 GeV m4 = 120 GeV a/p? = 60 GeV 

(101) 

7712 = 80 GeV ms = 200 GeV mg = 300 GeV ^Jpl = 20 GeV. 

For the decaying particle we vary the value of M = ^/p^ in a range between M = 100 GeV 
and M = 800 GeV and plot the real and imaginary part. The result is shown in Fig. |3 
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Figure 4: Real part (top) and and imaginary part (bottom) for the scalar planar two-loop 
three-point diagram of Eq. (|42|) with decay mass M between 100 GeV and 800 GeV. 
For the values of masses and momenta we use the standard set given in Eq. (jl01|) . 
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A three-particle threshold is expected for M = 1712 +^3 + 7/14 = 300 GeV. At this point 
we see that the imaginary part starts to differ from zero. The two-particle threshold at 
M = rrii + ms = 320 GeV is characterized by a sharp peak in the real part, accompanied 
by an instant decrease of the imaginary part, leading to a vertical tangent to the curve at 
this point. For the second two-particle threshold in this energy range at M = mi + m2 = 
500 GeV, the situation is the opposite. The imaginary part shows a sharp peak while the 
real part increases with vertical slope. 



5.3 The topology of the subtraction terms 

If we consider tensor integrals we have to take subtraction terms into account. As it was 
shown in the previous section, for the planar two-loop three-point master diagram with 
power {ko — ki)"^ we need a single subtraction. Each subtraction replaces the propagator 
factors by the subtracted ones. In the present case one has 



{ko - kifd'^kdH / P1P2 \ _ f{ko- kiYdVl r {ko - kifd^kdH 
1\IVW\P^ V ~ Pi,iP2,i) ~ J iww^K ~ J l\J\Jw^Pe 

(102) 

(for higher order subtractions the scaled propagator factors appear in higher powers). If we 
perform the limit k, —>■ for Pi 1 = {k + npiY — ml i—irj and P2,i = {k — ^^2)^ — ^2 1 ^ 
the propagator factors loose their dependence on the outer momenta Pi. In this case the 
momentum scheme given for the original diagram is no longer valid. We first have to 
perform a partial fraction decomposition for the subtraction term. 



{ko- kiYd^kd^l _ ] 
Pi,iP2,iPzPaP^P% ~ rn\i - m 



1,1 



(fco - kifd^k d^l f {ko - kifd^k d^l 



Pi^iP^PaP^Pq J P2,lP'iPAPbP& 

(103) 

After that, the lines can be rearranged according to their momenta. The planar diagram 
is shown in Fig. El together with the subtraction. For convenience, relative factors are not 
mentioned in the diagrammatic representation. 



6 Numerical integration 

After having performed all other integrations analytically, only the two integrations over 
/cq and Vq are left. The integration region is given by the non- vanishing of the factors S^'^ 



30 



p 




n=l i=l 



EE 



2 jk 



P 




Figure 5: Planar topology T° (left hand side) and the subtraction terms (right hand side). 
For convenience, relative factors are not mentioned in the diagrammatic representation. 
The numerator factor for the master integral and the subtraction terms read (/eq — ki)". 



K[m = -4vr^ [^(±(^0 + iom-en)ei-c^m) + e{±{K + iomOn)ei<Pm)] (104) 



depending on the linearizations fco = feg ± fci and Iq = I'q ± li. If we take into account that 



the integration regions can be seen to be triangles which are bound by the off-diagonal 
kg + I'q = and the different unequalities coming from 6n and (pm- In case of the index 
couples {n,m) = (1, 4), (2, 5), (0, 6) the integration region vanishes identically. The non- 
vanishing integration regions are shown in Fig. IHl in case of the linearizations fco = + 
and Zo = ^0 + ^1- Two remarks are in order at this point. 

• Special care has to be taken for the case (n, m) = (0, 6). In this case, the integration 
region vanishes only in the limit /tj, Aj —>■ 0. However, we have thoroughly checked 
that in this limit the integrand does not diverge, i.e. the integral does not give a 
non- vanishing contribution even with a vanishing integration region [U^ . 

• In case that one of the masses vanishes, integration regions might vanish. In Fig.lHlwe 

see that if Ei = q^, i.e. if pf = 0, the triangles (n, m) = (1, 6), (0, 4) vanish. In case of 
the opposite linearization and with = 0, the same holds for (n, m) = (2, 6), (0, 5) 
instead. However, the integrand does not diverge in these cases. 



which read 



04 = ±2(/[, -Ei± q,), 05 = ±2{1'q + E2± g.). 




(105) 
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Figure 6: The integration regions for the last (numerical) two-dimensional integration are 
triangles which are bound by the line /cq + Zq = 0, a vertical and a horizontal line, according 
to the condition that in Eqs. (IUH|) or (jHH) does not vanish. The integration regions 
for different values of index couples (n, m) are shown for the linearizations = k'^ + ki 
and Iq = Iq + li. In case of the opposite linearizations, the signs in front of the change. 
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The integrations are performed using numerical routines like the Monte Carlo integration 
routine VEGAS [13131] and the program Divonne of the library CUBA jlH]. While nu- 
merical routines sometimes overestimate the precision of the result IM] , the use of two 
different routines enables us to have a measure for the accuracy of the result. In this section 
we will present a few examples in order to demonstrate the reliability of our method. 

6.1 The decay via top quark loop 

As a first example we consider the basic integral in Eq. (jl(J2p occuring in a process where 
the boson couples to a triangle with top quark current (cf. Fig. [7j). We have started 
the discussion already in Eq. ()102|) and will continue at this point. The masses are chosen 
to be 

mi = m2 = m3 = m^, = rrir, = mz, and mg = mi (106) 

where we use the values rrit = 178 GeV, rriz = 91.1876 GeV and mi = 0. The 
subtraction masses which should have different values are taken as mi^i = 150 GeV and 
7712,1 = 160 GeV. Because the squared mass of the Z° boson is below all possible thresholds, 
the imaginary part of this integral vanishes. Even though physical masses are chosen in 
this paragraph, we stress that instead of calculating the whole process we are dealing with 
the occuring planar basic diagram only. 

We calculate the integral for 50 points between ki = 0.01 and ni = 0.5 for the (di- 
mensionless) parameter ki. As a next step we approximate the points by a low degree 
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degree 


real part 






0"Rc 






Xlj{n - 1) 





-2.941141 


10- 


-4 


1.4 


10- 


-9 


1.2 ■ 10^ 


1 


-3.484637 


10- 


-4 


2.6 


10- 


-9 


3.8-105 


2 


-3.361064 


10- 


~4 


3.9 


10- 


-9 


1.4- 10^ 


3 


-3.361950 


10- 


-4 


5.2 


10- 


-9 


1.4-10° 


4 


-3.361836 


10- 


-4 


6.8 


10- 


-9 


1.3 - 10° 


5 


-3.361884 


10- 


-4 


8.6 


10- 


-9 


1.3 - 10° 


6 


-3.361693 


10- 


-4 


1.1 


10- 


-8 


1.1 - 10° 


X 


-3.365555 


10- 


-4 


7.2 


10- 


-9 





Table 1: Extrapolation of the contributions for V2 using VEGAS 



degree 


real part 






0"Rc 






Xlj{n - 1) 





-2.81249 


10- 


-4 


3.9 


10- 


-8 


1.7- 10^ 


1 


-3.52842 


10- 


-4 


8.7 


10- 


-8 


5.0- 102 


2 


-3.36111 


10- 


'A 


1.3 


10- 


-7 


1.0-10-1 


3 


-3.36265 


10- 


-4 


1.9 


10- 


-7 


7.4-10-2 


4 


-3.36301 


10- 


-4 


2.5 


10- 


-7 


7.3 - 10-2 


5 


-3.36294 


10- 


-4 


3.2 


10- 


-7 


7.3 - 10-2 


6 


-3.36209 


10- 


-4 


4.1 


10- 


-7 


7.1 - 10-2 


X 


-3.36530 


10- 


-4 


3.3 


10- 


-7 





Table 2: Extrapolation of the contributions for V2 using Divonne 



polynomial in order to extrapolate to ki = 0. The result of this extrapolation can then be 
compared with the result for the integral where we used Hi = from the very beginning. 
For the integration with VEGAS, 260 000 internal points were used, whereas for Divonne 
we use a standard precision of 10-^. The results of this fit are shown in Tables Q and |21 If 
Ki = is used from the very beginning and all integrations are done analytically except for 
the last two ones, we obtain the results shown in the last lines of Tables ^ and El marked 
by a "x". Obviously, these values are already reliably approximated by the second order 
fit. In case of VEGAS, the quantity /{n — 1) for n sampling points as a measure for 
the reliability of the polynomial fit of the numerical integrations approaches the optimal 
value 1.0 quite fast. The deviation of the order of 1 in case of Divonne is not problematic 
because we use this second method only to check the VEGAS calculation. The result is 
shown graphically in Fig. |S1 In this figure we show the results of the numerical calculation 
as a function of ki together with the result for the semi-analytical calculation for ki = 0. 
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0.1 0.2 0.3 0.4 0.5 



Figure 8: Real part of the integral V2 in Eq. ()102|) as a function of the parameter Ki. 
For the masses and momenta we used values related to the decay of the boson into 
a top-quark loop, as explained in the text (cf. Fig. [7j). Shown are the results for the 
numerical calculation (curve) and the result for the semi-analytical calculation with ki = 
(horizontal straight line). In both cases we have used VEGAS. 
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degree 


real part ctrc Xnc/i^-'^) 


imag. part aj.^ xL/i^-'^) 


2 
3 


2.71015-103 3.0-10-1 8.5-10-1 
2.71047-103 4.2-10-1 8.3-10-1 


6.67079-103 7.7-10-1 1.1-10° 
6.6712-103 1.1-10° 1.1-10° 


X 


2.71125-103 7.0-10-1 


6.6703-103 1.7-10° 


Table 3: Extrapolation of the contributions for V21 using VEGAS 


degree 


real part (Trc Xrc/(^-1) 


imag. part aim Xim/ ~ 


2 


2.7094-103 1.0-10° 4.4-10-1 


6.6691-103 2.7-10° 3.7-10-1 


X 


2.7042 - 103 2.6 - 10° 


6.6736-103 6.6-10° 



Table 4: Extrapolation of the contributions for V21 using Divonne 



The error bars for the latter are plotted at different points in order to allow for a better 
comparison. Also optically the results for numerical and semi-analytical calculation match 
very well. This gives us confidence in the correct implementation of the algorithm. 

6.2 Rotated reduced planar topology 

Next we consider the numerical results for a specific example of the rotated reduced planar 
topology T^. The physical starting point is a decay of the Z° boson via a lepton loop, 
coupling to the outer lepton and anti-lepton legs via W bosons. The integral 

contains two subtraction terms, for the fc-loop as well as for the /-loop. In this case we 
consider the limit ni = 0.8 - Ai — > 0. The diagram with rotated reduced planar topology is 
shown in Fig. together with the three different types of subtraction terms. The masses 
take the values mi = 1712 = = mi = 0, = m^^ = 80.425 GeV IHSj, and mg = mi = 0. 
The momenta are given by ypf = -^pf = mi = and -y/p^ = mz- For the three subtraction 
masses we take mi 1 = 200 GeV, m2,i = 300 GeV, and m4^i = 110 GeV. In Table 121 we show 
the results for the numerical integration with VEGAS for the real and the imaginary part. 
Only the values at polynomial degrees 2 and 3 are shown. They are in good agreement with 
the result at polynomial degree 2 obtained by using Divonne (Table |3)). The results for the 
semi-analytical calculation are again given in the last lines of the two tables. Looking at 
the values obtained by using VEGAS in Fig. [HI one sees good agreement. 
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0.1 0.2 0.3 0.4 0.5 

Ai 




0.1 0.2 0.3 0.4 0.5 



Ai 

Figure 9: Real part (top) and imaginary part (bottom) of the integral V|i in Eq. ()1()7|1 as 
a function of the parameters Ai and ki = 0.8 ■ Ai. Shown are the results for the numerical 
calculation (curve) and for the semi-analytical calculation with ki = Ai = (horizontal 
straight line), both obtained by using VEGAS. 
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Figure 10: Rotated reduced planar topology (first diagram) and the subtraction terms. 
For convenience, relative factors are not mentioned in the diagrammatic representation. 
Summations over n = 1,2, i = 1, . . . , jk and j = 1, . . . ,ji are assumed (cf. Fig. Ej). The 
numerator factor for the master integral and the subtraction terms reads (A;o + fci)"(^o + ^i)^- 

6.3 Results for varying decay mass 

For the integral with two subtractions for the fc-loop, 



we vary the decay mass M = \fp^ in order to identify different thresholds. The values for 
physical masses and outgoing momenta are again taken from Eq. (jlOlj) . for the subtraction 
masses we use mi i = 100 GeV, m2,i = 200 GeV, mi 2 = 350 GeV, and m2,2 = 450 GeV. In 
varying M we can observe the behaviour of the real and the imaginary part. The behaviour 
is shown in Fig. ^2 for decay mass values close to M = 320 GeV. As in case of the simple 
example shown in Fig. |3] we notice a sharp peak for the real part and a vertical slope 
for the imaginary part at the point M = 320 GeV which corresponds to the two-particle 
threshold = = {1714 + ms)^. In a second step we are looking more closely at the 
threshold region. Taking the value M = 325 GeV, we analyze the real and imaginary parts 
for K2 = 0.6 ■ Ki — > 0. In Table we show results of the polynomial fit for the results 
obtained by VEGAS. The last line of Table El and Table IHl are the semi-analytic results 
obtained by using VEGAS and Divonne, respectively. The results are shown in Fig. ^1 
Again both methods lead to the same result in the limit K2 = 0.6 ■ ki — 0. 




(108) 
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Figure 11: Real part (top) and imaginary part (bottom) of the integral V3 in Eq. ()108|) as 
a function of the decay mass M = i/p^ close to M = 320 GeV. The values for the masses 
and outgoing momenta are taken from Eq. ()101|) . 
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Figure 12: Real part (top) and imaginary part (bottom) of the integral V3 in Eq. ()1()8|1 as 
a function of the parameters ki and k,2 = 0.6 ■ ni for the decay mass value M = a/p^ = 
325 GeV. Shown are the results for the numerical calculation (curve) and for the semi- 
analytical calculation with /t2 = ki = (horizontal straight line), both obtained by using 
VEGAS. 
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degree 


real part (Jrc xL/(^-1) 


imag. part oi,,, xlJ {n - 1) 


3 
4 


-3.00814-10-1 1.2-10-5 4.2-10° 
-3.00692-10-1 1.5-10-5 8.9 - 10"! 


-1.519135-10-1 5.2-10-6 3.6-10° 
-1.518660-10-1 6.8-10-6 1.2-10° 


X 


-3.00689-10-1 1.7-10-5 


-1.51865-10-1 1.4-10-5 


Table 5: Extrapolation of the contributions for V3 using VEGAS 


degree 


real part ctrc xL/ " 1) 


imag. part aim Xim/ - 1) 


X 


-3.0083 - 10-1 2.9 - 10"^ 


-1.5196-10-1 1.4-10-^ 



Table 6: Extrapolation of the contributions for V3 using Divonne 



6.4 Combined integration 

After having performed these tests, the semi-analytical evaluation of the partial results from 
the different master integrals appears to be feasible. This might be done by performing the 
numerical integrations separately for each master integral and afterwards combining the 
results. However, it turns out that the more stable method is to first perform the analytical 
integrations separately and then to do the numerical integration for the sum of all basic 
integrands relevant for a given process. This is especially true for squared momenta on the 
threshold or close to it. As an example we use 

After the tensor reduction explained in Sec. 3, we are left with the original planar topology 
and the two rotated reduced planar topologies, both with different powers of (/cq T ^1) and 
(^0 -F^i); and diagrams of simpler topologies. The master integrals for the non-reduced and 
rotated reduced planar topologies can be calculated by applying the subtractions explained 
in Sec. 4. For the convergent part the limit Kj, Aj — can be performed. All integrations 
up to the last two can be done analytically. The combination of the different UV-finite 
contributions with planar and rotated reduced planar topology at this point is called Vex- 

As explained before, there are two possibilities to proceed. We can combine all ana- 
lytical results to a single expression and integrate this expression numerically (combined 
integration), or we can perform the numerical integration of the analytical expression for 
each of the master integrals separately and add up the results afterwards (separated inte- 
gration). We will perform both calculations in the following to compare them. 

For the physical masses and outgoing momenta we again use the values in Eq. ()101|] . 
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The square of the decay mass M = a/p^ = 500 GeV is chosen to be equal to the two-particle 
threshold = {m^ + mg)^. For the subtraction masses we take values with an offset mo 
which in this example runs from GeV to 200 GeV,^ 

mi 1 = 5 GeV + mn mi 9 = 25 GeV + mn 

(110) 

m2,i = 15 GeV + mo m2,2 = 35 GeV + mo- 

Subtraction terms and corresponding subtraction masses for the second loop momentum / 
are not necessary, since in this particular case the contributing planar and rotated reduced 
planar topologies do not contain sub divergences in /. The results are shown in Fig. El The 
left hand side of Fig. shows real and imaginary parts for the separated integration, the 
right hand side shows real and imaginary parts for the combined integration. It is obvious 
that the combined integration is much more stable. 

This result is easy to understand. Assuming that the relative variances are the same 
for all n partial integrals and for the complete integral, we can use the same number of 
intersection points in the Monte Carlo integration for each seperate integral and for the 
complete integral. In this case, the standard deviation of the seperately integrated and 
summed integrals is bigger by a factor of y/n than the standard deviation of the complete 
integral PjlMj. 



7 Conclusions 

Using the parallel/orthogonal space method, we have elaborated on a new method to 
calculate the planar and rotated reduced planar ladder two-loop three-point diagrams in 
the general massive case. Our method does not depend on a specific physical decay process 
but is suitable to calculate all planar two-loop vertex diagrams which may arise in the 
Standard Model of electroweak interactions using Feynman gauge. 

We have presented a tensor reduction which reduces the mentioned topologies with 
arbitrary tensor structure to a set of master integrals with strongly restricted, well-defined 
numerator structure. We have developed and implemented an algorithm for the semi- 
analytical calculation of the UV-finite part of those master integrals which still have the 
mentioned topologies after tensor reduction. While six of the eight integrations are done 
analytically, for the remaining two integrations we used different numerical methods. 

^Note that Vex still does not need to be independent of the subtraction masses because only in the sum 
with the finite parts of the subtraction terms this dependence will vanish. 
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50 100 150 200 50 100 150 200 

mo[GeV] mo[GeV] 

Figure 13: Real part (top) and imaginary part (bottom) of the integral Vex as a function 
of the mass parameter itlq of the subtraction masses for the decay mass value on the two- 
particle threshold at M = a/p^ = 500 GeV. The integral Vex is the sum of the UV-finite 
parts of the master integrals corresponding to %x in Eq. which after tensor reduction 
keep the planar and rotated reduced planar topology. Shown are the results for separated 
(left hand side) and combined numerical integration (right hand side) with VEGAS. 
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The integration of the UV-finite part of the important master integrals was checked 
exhaustively in the limit of vanishing subtraction parameters by comparing the semi- 
analytical calculation with the results of a numerical integration. As the final step of 
the semi-analytical calculation, the numerical integration turns out to be stable, i.e. one 
obtains the same result as for the above mentioned limit. This is valid especially in the 
case where the numerical integration is not done for each master integral separately but 
for the sum of all analytical results necessary for the given process. 

For the remaining master integrals and subtraction terms with simpler topologies there 
exist other methods. In part, these methods are already implemented in ;^[oops. Therefore, 
in writing links to existing codes, these cases can easily be implemented into our calcu- 
lations. The still missing topologies are the two-loop three-point topologies containing a 
two-point subloop. These parts still have to be implemented |24| I25j . After having done 
this, we will be able to evaluate semi-analytically all IR-finite NNLO corrections containing 
planar two-loop vertex topologies as they arise in the Standard Model. 

Finally, because the calculation of the scalar integrals for the non-planar two-loop vertex 
function 1^ is quite similar to the calculation of the scalar integrals for the planar two-loop 
vertex function j^, our algorithm is expected to be easily extendable to the non-planar 
two-loop vertex topology. Having extended our algorithm in this direction, we are then 
able to calculate all NNLO vertex corrections as they arise in the Standard Model. 

As far as possible, the computer codes written in the course of this paper are developed 
and implemented for most general decay processes. Because of this, we plan to include our 
calculations in the ^^Coops-pioiect in the near future. 
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A Integral basis 



In this appendix we list the analytical results for the basic integrals in Eq. in terms 
of the parameters r^, rf, tq, sq, and to. In order to calculate these integrals, we introduce 
cutoffs and Aj. The results depending on these cutoffs are extremely long and contain 
products of two logarithms as well as dilogarithms [HI]- However, one can show that these 
divergences cancel in the sums in Eqs. fl()7|l and (jHH|) . In order to make the reader familiar 
with these features, we explain the calculation for the simplest example. 



A.l The integral JFqq 



The regularized integral is given by 



0,0 



dt 



ds 



1 



(ro — TsS — Ttt + iriY — 4:St 



TsTo + (2 - rsrt)At - A^ - r^A 
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2(rsro + (1 - rsrt)At) 

TsTo 
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A, + 
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rtTo 



TtTo + (2 - rsrt)As - rfAt - nA 



-In |ro(l + rsn) + (1 - rsrt){rsAs + nAt + A)| - In 



2{rtro 
+ ivr^(ro) 
'"o 



rsrt)As) 



2(1 - r,rt) 



where 



A := ^ {ro - TsA, - nAtf - AAsAf 



(Al) 



(A2) 



What to do with the singularity? As a first step we can separate the divergent part which 
depends on A^, At and A from the convergent part. But how does the divergent part cancel? 
If we look at Eq. (j88|) where JFJ q appears, we see that it actually appears as part of a sum. 
If we concentrate on the summation over n and consider how Mii ^ and A/Jj ^ are defined 



6,j - iij - -Pi,j\ 



in Eqs. fj53|) and ()54|) . and if we take into account that P2, 

Therefore, an expression which is of relevance for the subtracted integrals is the sum 



we see that U^^^^ = -Af^i^^- The same is valid for TV^^^-b and A^^^^. ;,, i.e. A/].^.;, - -^^ej^b- 



E E 

n=l m=xJ 



\n+in rst 



•^0,o(^S' "^ti '^ri,i,m,j)- 



(A3) 
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In replacing by A°/eA and At by A°/eA with ca ^ 1 with arbitrary but positive and 
fixed values A° and A? we indeed can show that 



E E (-l)"+'"-^oT^(r.,n,r„,)=0(eA) 



(A4) 



^ ^\ n+m ^st ,div / 

n=l m=xfi 

For eA 0, therefore, the divergent part will cancel in the difference. We only have to 
consider the convergent part 



■J f 



st,conv 
0,0 



'^O 



TsTt 



In 



'^o 



2 1 - r^rt) 



Finally, we see that 



E E (-1) 



n+m 



^ E (ei. - 6.) = 



(A5) 



(A6) 



n=l m=xfi 

(cf. Eq. (jyij) including remarks given there). Therefore, we can skip all terms which are 
linear in tq. The result which we have to implement is 



^sijimpl 
0,0 



'^O 



[iTc9{ro) - In |ro|) . 



(A7) 



1 - TsTt 

In the following we will only give results which are to be implemented. For the integrals 
with one index we will show only results for JF*. Results for JF^ can easily be obtained by 
replacing sq ^ to and rt. 

The integral JFq q is used for illustrative reasons only. It does not appear in calculations 
because the integration region vanishes. The same is valid for the integrals J-'fj. Therefore, 
it turns out that the integrals JF^ and J^^ are sufficient for calculations within the Standard 
Model of electroweak interactions using Feynman gauge. 



A. 2 The integrals JT: 
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+ ((2-r,n)to + r,ro)^(to)ln 
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(A9) 



5 -r-«,impl 
I s-' 2 



{1 / r^ro \^ , 2 2\ 

((6 - Gr^rt + rlrt)to + 2(3 - rsTtYsTf^ | In 



1 - TsTt 



ro\ 



6 - Qrsn + r^r^ )to + 2(3 - rsTtYsrotQ + r^rQ j Re Li 1 



to(l - r^rt) 



-«vr<^^(ro) 



2 V 1 - rsTt 



(6 — Gr^r^ + r^r^^ 



(6 - Qrsn + r^r^ )to + 2(3 - TsTtYsTQ 



1 - r^rt 

+ ((6 - Qr^rt + r^ri2)to + 2(3 - rsrt)rsroto + r^r^) In 



1 + 



rsro 



-e{to) ((6 - 6r,rt + r.V^)^^ + 2(3 - r,rt)r,roto + r^rl) In 



to(l - r^rj) 
rsro 



to + 



1 - r^n 



(AlO) 



References 

[1] W. Hollik and G. Duckeck, "Electroweak Precision Tests at LEP", Springer 2000 

[2] D.Y. Bardin, J. Phys. G29 (2003) 75 

[3] W. Hollik et ai, Acta Phys. Polon. B35 (2004) 2533 

[4] [LEP Collaborations] , "Precision electroweak measurements on the Z resonance" , 
Report No. CERN-PH-EP-2005-041, ,arXiv:hep-ex/0 509008|, 
submitted to Physics Reports 

[5] J. A. Aguilar-Saavedra et al. [ECFA/DESY LC Physics Working Group], 

"TESLA Technical Design Report Part IIL Physics at an e^e~ Linear Collider" , 
Report No. DESY-TESLA-2001-23 farXiv:hep-ph/0106315 



[6] S. Heinemeyer and G. Weiglein, "Electroweak precision tests with GigaZ", 

published in Batavia 2000, Physics and experiments with future linear e'^e~ colliders, 
pp. 511-514, Report No. LC-TH-2001-001 | arXiv:hep-ph/0012364| 



47 



[7; 



[9 
[10 

[11 

[12 
[13 
[14 

[is: 



G. Weiglein, Eur. Phys. J. C33 (2004) S630 

G. Weiglein et al. [LHC/LC Study Group], "Physics interplay of the LHC and the 
ILC", Report No. SLAC-PUB-10764 |arXiv:hep-ph/0410364^ 

S. Heinemeyer, "Loop calculations: Summary", 

Report no. CERN-PH-TH-2004-159 [arXiv:hep-ph/0408269l 

D. Kreimer, Phys. Lett. B292 (1992) 341 

N.I. Usyukina and A.I. Davydychev, Phys. Lett. B348 (1995) 503 

A.I. Davydychev and J.B. Tausk, Nucl. Phys. B465 (1996) 507 

J. Fujimoto, Y. Shimizu, K. Kato and T. Kaneko, Int. J. Mod. Phys. C6 (1995) 525 

A.A. Pivovarov, Phys. Atom. Nucl. 63 (2000) 1646 [Yad. Fiz. 63N9 (2000) 1734] 

K.G. Chetyrkin, R. Harlander, J.H. Kiihn and M. Steinhauser, 
Nucl. lustrum. Meth. A389 (1997) 354 



[16] S. Groote, J.G. Korner and A.A. Pivovarov, Phys. Lett. B443 (1998) 269; Nucl. Phys. 
B542 (1999) 515; Eur. Phys. J. Cll (1999) 279; Phys. Rev. D60 (1999) 061701; 
Phys. Rev. D61 (2000) 071501; Phys. Rev. D65 (2002) 036001; Eur. Phys. J. 
C36 (2004) 471; "On the evaluation of a certain class of Feynman diagrams in 
x-space: Sunrise-type topologies at any loop order". Report No. MZ-TII-05-08 
|arXiv:hep-p h/0506286|; S. Groote and A.A. Pivovarov, Nucl. Phys. B580 (2000) 
459 

[17] A. Ghinculov and Y.P. Yao, 

Nucl. Phys. B516 (1998) 385; Phys. Rev. D63 (2001) 054510 

[18] A.I. Davydychev and V.A. Smirnov, Nucl. lustrum. Meth. A502 (2003) 621 

[19] U. Aglietti and R. Bonciani, Nucl. Phys. B668 (2003) 3; B698 (2004) 277 

[20] W. Bernreuther, R. Bonciani, T. Gehrmann, R. Heinesch, T. Leineweber, P. Mastrolia 
and E. Remiddi, Nucl. Phys. B706 (2005) 245; B712 (2005) 229 



48 



[21] R. Bonciani, P. Mastrolia and E. Remiddi, Nucl. Phys. B661 (2003) 289 
[Erratum B702 (2004) 359]; B676 (2004) 399; B690 (2004) 138; 
R. Bonciani, A. Ferroglia, P. Mastrolia, E. Remiddi and J.J. van der Bij, 
Nucl. Phys. B701 (2004) 121 



[22 
[23 
[24 



[25 
[26 
[27; 
[28 
[29 

[30 

[31 



[32 

[33 
[34; 
[35 



M. Caffo, H. Czyz and E. Remiddi, Nucl. Phys. Proc. Suppl. 116 (2003) 422 

P. Mastrolia and E. Remiddi, Nucl. Phys. Proc. Suppl. 116 (2003) 412 

S. Bauberger, F.A. Berends, M. Bohm and M. Buza, Nucl. Phys. B434 (1995) 383; 
S. Bauberger, M. Bohm, G. Weiglein, F.A. Berends and M. Buza, 
Nucl. Phys. Proc. Suppl. 37B (1994) 95 

P. Post, Ph.D. thesis, Mainz University, 1997 

P. Post and J.B. Tausk, Mod. Phys. Lett. All (1996) 2115 

H.S. Do, Ph.D. thesis, Mainz University, 2003 

G. Passarino, Nucl. Phys. B619 (2001) 257 

A. Ferroglia, M. Passera, G. Passarino and S. Uccirati, 
Nucl. Phys. B680 (2004) 199 

S. Actis, A. Ferroglia, G. Passarino, M. Passera and S. Uccirati, 
Nucl. Phys. B703 (2004) 3 

W. Hollik, U. Meier and S. Uccirati, "The effective electroweak mixing angle sin^(^eff) 

with two-loop fermionic contributions". 

Report No. MPP-2005-77 arXiv:hep-ph/0507158| 

M. Awramik, M. Czakon, A. Freitas and G. Weiglein, 

Phys. Rev. Lett. 93 (2004) 201805; Nucl. Phys. Proc. Suppl. 135 (2004) 119 
D. Kreimer, Phys. Atom. Nucl. 56 (1993) 1546 

A. Czarnecki, U. Kihan and D. Kreimer, Nucl. Phys. B433 (1995) 259 
A. Frink, U. Kilian and D. Kreimer, Nucl. Phys. B488 (1997) 426 



49 



[36] D. Kreimer, Mod. Phys. Lett. A9 (1994) 1105 

[37] A.V. Kotikov, Phys. Lett. B259 (1991) 314 

[38] J. Fleischer, A.V. Kotikov and O.L. Veretin, 

Phys. Lett. B417 (1998) 163; Nucl. Phys. B547 (1999) 343 

[39] V.A. Smirnov and E.R. Rakhmetov, 

Theor. Math. Phys. 120 (1999) 870 [Teor. Mat. Fiz. 120 (1999) 64] 

[40] T.G. Birthwright, E.W.N. Glover and P. Marquard, JHEP 0409 (2004) 042 

[41] J. Fleischer, V.A. Smirnov, A. Frink, J.G. Korner, D. Kreimer, K. Schilcher and 
J.B. Tausk, Eur. Phys. J. C2 (1998) 747 

[42] A. Frink, Ph.D. thesis, Mainz University, 2000 

[43] G.P. Lepage, "Vegas: An Adaptive Multidimensional Integration Program", 
Report No. CLNS-80/447 

[44] T. Ohl, Comput. Phys. Commun. 120 (1999) 13 

[45] T. Hahn, "CUBA: A library for multidimensional numerical integration". 



Report No. MPP-2004-40 |arXiv:hep-ph/0404043| 
[46] D. Kreimer, Nucl. lustrum. Meth. A389 (1997) 323 
[47] L. Briicher, Nucl. lustrum. Meth. A389 (1997) 327 
[48] J. Franzkowski, Nucl. lustrum. Meth. A389 (1997) 333 

[49] A. Frink, J.G. Korner and J.B. Tausk, "Massive two-loop integrals and Higgs physics" , 
published in the Proceedings of the Joint Particle Physics Meeting, Ouranoupolis, 
Greece, May 1997, pp. 175-200 

[50] C. Bauer and H.S. Do, Comput. Phys. Commun. 144 (2002) 154 

[51] C. Bauer, A. Frink and R. Kreckel, Report No. MZ-TH/00-17 [arXiv:cs.sc/0004015l 

[52] D. Kreimer, Z. Phys. C54 (1992) 667 



50 



D. Kreimer, Int. J. Mod. Phys. A8 (1993) 1797 

A. Frink, Diploma thesis, Mainz University, 1996 

G. Passarino and M.J.G. Veltman, Nucl. Phys. B160 (1979) 151 

A. Andonov et al, "SANCscope - v.0.41", |a^iv:hep-ph/0411186| 



R. Harlander and M. Steinhauser, Prog. Part. Nucl. Phys. 43 (1999) 167 
U. Kilian, Ph.D. thesis, Mainz University, 1996 
L.D. Landau, Nucl. Phys. 13 (1959) 181 
M. Bohm, A. Denner, H. Joos, 

"Gauge Theory of the Strong and Electroweak Interaction", Teubner, 2001 
R.E. Cutkosky, J. Math. Phys. 1 (1960) 429 
M.M. Knodel, Ph.D. thesis, Mainz University, 2005 
R.H. Kleiss and A. Lazopoulos, 



"Error in Monte Carlo, quasi-error in quasi-Monte Carlo", arXiv:hep-ph / 0504085^ 



A.F.W. van Hameren, Ph.D. thesis, 
Nijmegen University, 2001 |arXiv:hep-ph/0101094 



S. Eidelman et al. [Particle Data Group], Phys. Lett. B592 (2004) 1 



51 



